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Abstract 

We treat a constant coefficient hyperbolic system in one space variable, 
with zero initial data. Dissipative boundary conditions are imposed at the 
two points x = ±1. This problem is discretized either by a spectral or 
pseudospectral approximation in space. We demonstrate sufficient conditions 
under which the spectral numerical solution is stable; moreover, these 
conditions have to be checked only for scalar equations. The stability 
theorems take the form of explicit bounds for the norm of the solution in 
terms of the boundary data. The dependence of these bounds on N, the number 
of points in the domain (or equivalently the degree of the polynomials 
involved), is investigated for a class of standard spectral methods, including 
Chebyshev and Legendre collocations. 
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INTRODUCTION 


In this paper we study the stability of spectral polynomial methods for 
the approximation of initial-boundary value hyperbolic systems with constant 
coefficients and dissipative boundary conditions (Assumption I). We show that 
any particular spectral method is stable when applied to a system, if it 
satisfies certain conditions (Assumptions II and III) for the corresponding 
scalar problem. 

Our treatment follows closely the approach of Kreiss [4], leading to an 
algebraic condition. However, in the spectral method the requirement that the 
approximate solution be a polynomial produces a different algebraic problem 
from the standard finite difference discretizations. 

We should also review briefly the stability results for spectral methods 
and hyperbolic systems. The stability of Chebyshev approximation of the 
scalar problem 



(u(x=l,t) = 0 

was proved in [2] for the Galerkin and tau methods, and in [7] for the 
collocation method. The same proofs hold for the Legendre method; they are 
based on energy estimates, showing that the norm 

1 2 
/ (1+x) w(x) u (x,t)dx 

-1 

is bounded by the data. Here the weight w(x) depends on the discretization 
used and equals (1 - x 2 )~ ^2 for Chebyshev, w(x) = 1 for Legendre. 

Reyna [8] suggested a stable numerical method for strictly hyperbolic 
systems with maximal dissipative boundary conditions, based on high mode 
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smoothing. His result applies to problems with variable coefficients. While 
his procedure ensures Z stability, our constant coefficients proof requires 
no smoothing; instead stability is obtained here in certain weighted 
Z 2 norms. 

The first results about hyperbolic systems with general dissipative 
boundary conditions were reported in [5], where it is shown that the numerical 
solutions decay in time, thus mimicking the differential solutions. 

In the present paper we combine the above results into a general stability 
theory for the stability of spectral methods. A future paper will discuss the 
convergence of the spectral approximations to the exact solution. 


Definitions 

Consider the first order hyperbolic system of partial differential 
equations 

3u 3 u 

= A , -1 < x < 1, t > 0. 

Here, 

u = u(x,t) = (u (1) (x,t),-«* ,u (n) (x,t))' 

is the vector of unknowns and A is a fixed nxn coefficient matrix. Since 
by hyperbolicity A is similar to a real diagonal matrix, we may assume 
without loss of generality that it is diagonal: 
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A = 


(1.1b) 


'A 

,0 




A 


I 




> 0 . 


The solution of this system is uniquely determined if we prescribe initial 
conditions 


(1.1c) 


u(x,0) =0, -1 < x < 1 


and boundary conditions 


(1 -Id) 


Jg I <— 1 * t) = L u U (-l,t) + g I (t) 
u II (l,t) = R u I (l,t) + g II (t), t > 0. 


In these formulas, 


s = s(t) ■ g II (t))' 


is a given n-vector, and 


(l.le) 


I r 
u = I u 

/V ** 


( 1 ) 


,U 


U)v 


r> 


11 r 

U = [u 


( 4 + 1 ) 


• • ,u 


(n)- 


is a partition of u into its inflow and outflow components — corresponding 
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to the partition of A in (1.1b) - while L and R are constant 

reflection matrices of order £ x (n-£) and (n-£) x £, respectively. 

The system (1.1a) - (l.le) is a well-posed problem in the sense that it 
satisfies an a priori energy estimate which we now describe. Define the 
spatial norm associated with a pair of positive weight functions 
u(x) = (oi I (x) ,m I][ (x)) , 

(1.2a) Hu(x)ll 2 = !iu(x)ll 2 = / lu^x)! 2 w^(x)dx + J |u^(x)| 2 m^*(x)dx. 

-1 -1 

Here and elsewhere in the paper we denote by |v| the Euclidean norms of a 
vector v; similarly |a| = sup|Av|/jv|. 

Then, for well-posedness one must obtain the following inequality, with some 
U 0 > 0: 

(1.2b) 

oo oo 

(q ~ n 0 ) / e -2(n " n 0 )t: llu(x,t)ll 2 dt < const / e -2 ^ n “ 1 V t |g(t) | 2 dt 

t=0 ~ t=0 

for all q > n Q . 

The parameter q^ measures the exponential time growth of the solution. We 
are interested in the case where there is no growth in time, and therefore 
postulate the following dissipativity requirement: 



Assumption I, There exists a constant 6 > 0 such that 


(1.3) 


R| • | L | < 1 - 6 < 1. 


The inequality (1.3) guarantees that waves originating at one of the 
boundaries are not amplified when reflected at the other one. Consequently, 
there is no time growth; in Appendix A we show that there is a weight to such 
that for any q > 0 

CO 00 

(1.4) q / e nt llu(x,t)ll 2 dt < const / e -2r)t |g(t)| 2 dt, 

t=0 ~ t=0 

i.e., tIq in formula (1.2b) may be taken as zero in this dissipative case. 

We study polynomial spectral and pseudospectral discretizations of 
(1.1). In any such approximation, one seeks a vector of N-degree polynomials 


such that 


v = v N (x,t) = (v^ 1) (x,t),...,v^ n) (x,t))' 


(1.5a) 


9v 

dt 


3v 

A + I* 


Here, Q(x) is a diagonal matrix of the form 


(1.5b) 


, q I (x) I 


Q(x) = 


£x£ 


„ 0 ) 

“ U) V«)x(. -1)/ 
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where q*(x), q^^Cx) are the N-degree polynomials which characterize the 
specific (pseudo)spectral method employed, and 


f V f I 11^, 

T = t(t) = (t , T ) 


is an n-vector to be determined by the set of boundary conditions (l.ld): 


(1.5c) 


! v I (— 1 , t > = L v n (-l,t) + g*(t) 
v II (l,t) = R v X ( 1 , t ) + g I][ (t). 


We shall call q*(x), q^(x) the forcing polynomials, since they appear as 
inhomogeneous terms in the discretization of the originally homogeneous system 
(1.1). Some examples are in order. 

In the spectral Galerkin Chebyshev method, one has: 


( 1 . 6 ) 


t i N 

q (x) -iT 0 (x) + l (-l) k T, (x) 

Z U k=l k 

II 1 ^ 

q (x) = j T q (x) + l T R (x) 


with T k (x) denoting the k-th Chebyshev polynomial. 

In the pseudospectral Chebyshev method, one may collocate only at the 
interior extrema of Tjj + }(x), yielding: 


q I (x) = q n (x) = T' +1 (x) 


(1.7) 
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or use as collocation points these extrema together with the downstream 
boundaries, in which case: 

f q X (x) - (l-x) T'(x) 

(1 * 8) II 

(q (x) = (1+x) T'(x). 

Similarly, in the pseudospectral Legendre case one has 
(1.9) q I (x) = q IX (x) = p ^ +1 (x) 


or alternatively 
( 1 . 10 ) 


/ q I (x) = (l-x) P'(x) 
( q IX (x) = (l+x) p'(x) 


with P^Xx) denoting the k-th Legendre polynomial. These and other examples 
are outlined in [2,3]. 

In this paper we provide a stability study for the (pseudo)spectral method 
(1.5). We define stability in terms of an a priori energy estimate analogous 
to the differential one: 


Definition (1.1) . (Stability). The approximation (1.5a) - (1.5c) is 
stable if there exist a weighting pair m(x) and constants a and > 0 

such that for all D > Dq we have 


( 1 . 11 ) 


“ -2(q-n n )t , 00 -2(n-n n )t 

(n“H n ) / e II v (x,t)ll dt < const N / e |g(t)| dt. 

t=0 t=0 
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We note that the exponential time growth factor here, tiq> need not be the 
same as the one in the differential estimate (1.2b). The other constant, a, 
accounts for a possible algebraically increasing dependence on the 

discretization parameter N. 

To present our results, let us introduce into (1,5) a new variable 
6 ~N^ x,t ^ anc *> following the procedure of [4], Fourier-transf orm the 

resulting equations with respect to time. Denote by £ the real dual 
variable corresponding to t, and by 


( 1 . 12 ) 


( i 5 v N (x,s) =^(e _rit v N (x,t)) 

( M = s(s) = 9 {e nt jg(t)) , S = n + i£ 


the transforms of e ^ v^ and e r|t g, respectively (these functions are set 
to zero for t < 0). The resulting equation for v = fv*(x,s), v^ 1 '(x,s))' 

cv '/v rw ' 

is then: 


(1.13a) 



I dv T 

= A di- + q (x) 


i 

T 


= A 


II 


dv 


II 


dx 


+ q Ii: (x) 


II 

T 


(V 


where the n-vector t = t(s) = is to be determined from the 

boundary conditions : 


( y X (-l,s) = L v n (-l,s) + ^(s) 
(v n (l,s) = R v I (l,s) + g«(s). 


(1.13b) 
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Using Parseval's relation we conclude: 

LEMMA 1.2. The approximation is stable if there exist a weighting pair 
oj(x) and constants a and ri^ > 0 such that for all s with 
Re s = n > n q the following holds 

(1.14) (n - ri 0 )llv N (x,s)ll 2 < const N 2a | jg( s ) | 2 . 

In sections 2 and 3 we discuss necessary conditions for the stability 
estimate (1.14). One such condition is the obvious requirement that the 

(pseudo)spectral method must be stable for a scalar problem. Another one is 
that the homogeneous two-point boundary value problem (1.13) should have no 
eigenvalues with Re s > iIq* In section 4 we show that a strengthened version 
of these necessary conditions guarantees stability. The novelty of our 
sufficient stability criterion is that it deals exclusively with the 
properties of the scalar model and there is no need to consider the 
complicated coupling of such scalar equations through the boundary conditions 
(1.5c). 

In the last sections we demonstrate various important cases in which our 
sufficient conditions for stability are satisfied and, thus, deduce stability 
theorems for pseudospectral discretizations of the initial-boundary value 
system (1.1). 
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2. THE SCALAR PROBLEM 

The system (1.13a) decouples into n scalar equations of the form 


(2.1a) 


sv 



q(x)i, 


-l < x < l . 


A A 

Here v = v N (x,s) and a, stand for any of the corresponding components of 

A 

Vn(x,s) and the diagonal of A, q(x) equals either q*(x) for a < 0, or 
q II (x) for a >0, and x = x(s) is to be determined by boundary conditions. 
In the scalar case, these boundary conditions amount to prescribing the 
upstream values: 


(2.1b) 


Vjq ( — 1 , s ) = g(s), a < 0, 

v N (l,s) = g(s), a > 0. 


If formula (1.14) is valid, there must be a similar stability estimate 
holds for (2.1). We make the necessary 

Assumption II . (Scalar stability). There exist a weighting pair w^(x), 
o) II (x) and constants a and q Q > 0, such that for all s with 
Re s = n > n Q , we have 

(2.2a) (q - q Q ) llv N (x,s)ll 2 I < const N 2a | v N (-l ,s) | 2 , a < 0, 

w 

(2.2b) (q - q Q ) |] v N (x, s ) II 2 IX < const N 2 “ |v N (l,s)| 2 , a > 0. 

0) 



3. THE EIGENVALUE PROBLEM 

Connected with (1.13) is the following eigenvalue problem: 


Let be the space of all polynomials of degree < N. We say that 

<{> £ 1% is an eigenfunction of (1.13) corresponding to an eigenvalue s, if 

& = = is a nontr l v ial solution of 


(3.1a) 



I d< ^N 1 . . I 
= A _ + q (x)t 


j.II 

_ II II. . 

- A + q (x) 


II 

T 


with boundary conditions 


(3.1b) 


^ I (-l,s) = L^ n (-l,s) 
. & I 3 : (l,s) = R^d.s) . 


In order to determine whether s * 0 is an eigenvalue, we proceed as 
follows. Equations (3.1a) form a system of ordinary differential equations 
in x, depending on a parameter s 


(3.2) 


( Sl nxn ‘ A * Q (x) I* 


This can be solved by a formal series 


£(x.) = ( sI " A |^) 1 = l s -k_1 A k Q (k) (x)x. 

k=0 


(3.3) 
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We note that since Q(x) consists of N-degree polynomials, the sum on the 
right include only the first N+l terms. Rewriting this in the usual 
partitioned form, we have 


LEMMA 3.1. Let^ P I (x,s) = p I (x,s;A I ) 
denote the diagonal matrices 


and P n (x,s) = p^(x, s ; A**) 


(3.Aa) 


i p 1 (x,s) 
P Ii; (x,s) 


N 

I 

k=0 


-k-1 


N 


I 

k=0 


-k-1 

s 


(q I (x))^ k) (A T ) k 
(q II (x))< k) (A 11 )". 


Then the polynomial solution of (2.2a) is given by 


(3.4b) 


^ I (x,s) = P I (x,s)t I 

II, . „II, . II 
V& (x,s) = P (x,s)x . 


Substituting this expression into the boundary conditions (3.1b), we end 
up with the homogeneous linear system: 


(3.5) 


P (-l,s) 


-RP i (l,s) 


-LP II (-1 ,s) 


P Ii: (l,s) 


= 0 . 


II 


Denoting the coefficient matrix by D(s) 
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(3.6a) 


D(s) 


( P I (-l,s) 
-RP I (l,s) 


-LP 


II 


.II 


(- 1 , 8 ) 

( 1 , 8 ) 


y 


we arrive at 


USMHA 3.2. A complex number s * 0 is an eigenvalue of (3.1), if and 
only if . 


Det [D(s)J = 0. 


Suppose now that an eigenvalue s exists, with Re s = n > n Q > 0. Then 
the corresponding eigenfunction is a nonvanishing solution to the homogeneous 
problem (1.13) with g = 0, in contradiction to (1.14). For future reference 
we state this as 

1*®MMA 3.3. A necessary condition for the stability estimate (1.14) to 
hold is that the eigenvalue problem (3.1) has no eigenvalues s with 
Re s > n Q > 0, i.e. . 


Re s = n > tIq > 0. 


(3.6b) 


Det[D(s)] t 0, 
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4. THE STABILITY ESTIMATE 

Here we present a stability theorem for hyperbolic systems with initial 
and boundary value conditions. 

In the previous two sections we exhibited two necessary conditions for the 
stability of (1.5). The first of these, (2.2), is a property of the 

individual scalar equations which constitute system (1.5a). The second one, 
(3.6), involves the coupling of these equations through the boundary 
conditions. Thus, in principle, in order to test this second condition, one 
should check the boundary detrminant of D(s) for each pair of the reflection 
matrices L and R. However, (3.6) is implied by a much simpler condition, 
which deals with the scalar equations separately and does not depend upon L 
and R. Moreover, this strengthened condition turns out to guarantee 
stability, as we shall presently prove. 

We make the following 

Assumption III . The polynomial solution of the scalar inflow problem: 

A I 

(4.1a) s v 1 = a + q*(x), a < 0 


satisfies 

(4.1b) |v I (l,s)| < |v I (-l,s)|, Re s = n > n 0 > 0. 


Similarly, in the outflow scalar case, 


s v 


'II dv 


II 


(4.2a) 


a 


dx 


, II, X 

+ q (x), 


a > 0, 
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we have 

(4.2b) |v n (-l,s)| < |v II (l,s)|, Re s = n > ri Q > 0. 

The estimates (4.1b) and (4.2b) are discrete analogues of differential 
estimates, as shown in Appendix A, cf. a. 11. 

We can now prove 

LEMMA. 4.1. Let Assumptions I, II, and III hold. Then, there exists a 
constant rig > 0 such that for all s with Re s = ti > iIq we have 

(4.3) Det[D(s) ] * 0. 


Proof . Consider the diagonal matrix P^(x,s) of (3.4a). Its j-th 
*(i) 

diagonal element p^ J '(x,s) is given by 


:(j> 


N 


( x , s ) = ( s “ a j 1 q I (x). 


That is, we have 


N 


^^x,s) 


(4.4a) 


P I (x,s) 


p^^x.s), 


^(i) A (1) 

where p J = p^ J (x,s) satisfies the inflow scalar equation 
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(4.5a) sp^ = - ' + q I (x), a j < °, 1 < j < £. 

Similarly, the diagonal matrix 

(4.4b) P IX (x, s) = 

consists of solutions to the outflow scalar equations 

* ( j ) 

(4.5b) sp (j) = + q II (x) a j > ^ < J < n. 

Fix s with Re s ~ n > tIq > 0. We first claim that both matrices P(-l,s) 
and P(l,s) are nonsingular. Indeed, if for some index j, 1 < j < 1, p^ 
vanishes on the left boundary 

p N j)(_1,s) “ °» 

then by the scalar stability Assumption I, see (2.2a), the solution of (4.5a) 
vanishes everywhere , 

ll^ ) (x,s)ll x = 0. 

U) 

This, however, contradicts formula (4.5a)# In the same way, one shows that 
P(l,s) is regular. Hence, the matrix D(s) can be factored (here the ratio 
of two diagonal matrices is understood as the diagonal of quotients) 
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(4.6) 


D(s) = 


r P (-1,8) 


,-RP (l,s) 


■LP II (-l,s)' 


P II (l,s) i 


P (-1,8) 



P (-1,8) 


P II: (l,s) 


Next we prove that the matrix E given by 


P T (l,s) 


-L 


(4.7a) E(s) = 

k 

P A (-1,8) 

is regular, with a uniformly bounded inverse: 


P Ii: (-l.s) 

P II (l,s) 


(4.7b) |E _1 (s)| < K = 1 + l R j— + , Re s = n > n 0 > 0. 


It is here that we make use of Assumptions I and III, Denoting 


r- r OJL§)_ 
P (-1,8) 


L = L 


P n (-l,s) 

P^l.s) 


(4.8) 
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we have, by (1.3), (4.1b) and (4.2b): 



Hence, E can be inverted 


(4.10a) 


E _1 (s) 




and the estimate (4.7b) follows: 



(4.10b) 


|e _ 1 (s) < (2 + Ir I 2 + III 2 ) 1/2 



Thus, the matrix D in (4.6) is regular, being the product of regular 
matrices, and the lemma is proved. 

The estimate (4.10b) is, in fact, the last ingredient needed for a 
stability proof for (1.5). We can state 


THEOREM 4*3. (Stability). Consider the hyperbolic system (1.1a) - 
(l.le), satisfying the dissipativity Assumption I. Then its (pseudo)spectral 
approximation (1.5a) - (1.5c) is stable, provided that Assumptions II and III 


are fulfilled. 
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Proof . System (1.5a) consists of l inflow equations 


(4.11a) 



d ;o> 

dx 


+ q I (x)x (j) 


1 < j < l 


and n-I outflow ones 


(4.11b) 



dv 


;<j> 


dx 




l < j < n. 


By the scalar stability Assumption II, each of these scalar equations is 
stable. That is, there exist a weighting pair m(x) = (w I (x) ,(o II (x ) )' and 
constants a and tIq > 0 such that for all s with Re s « q > we have 
for arbitrary vectors t 

(4.12a) (n - n 0 )«v^ ) (x,s)ll 2 I < const N 2ct |v^ j) (-l,s)| 2 1 < j < *, 

0) 

(4.12b) (n - q 0 )llv^ ) (x,s)ll 2 I] . < const N 2 “| v^ \l ,s) | 2 A < j < n. 

CO 

Using the spatial norm of (1.2a) 


llv N (x,s)ll 2 = llv N (x,s)ll 2 = llvj(x,s)l| 2 ] . + llvJ I (x,s)|| 2 I 

Oi 03 


II 


these inequalities can be added together and rewritten in concise form 


(n " V"~N (x ’ s)l, ‘ < COnst N 2 a [liN (-1 ’ s) | Z + 


'II, 


Res = n>riQ>0. 


(4.13) 
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It remains to estimate the two squared boundary terms inside the brackets on 
the right. To this end we rewrite our solution in the form, cf. (4.5) and 
(4.11) 


(4.14) 


vj(x,s) = P I (x,s)t I 
(x , s ) = P II (x,s)t 11 . 


Inserting this in the boundary conditions (1.5c) we find, as in (3.5), that 
t = t(s) is determined uniquely by 



I 1 

T 

Ai 

= D *(s) 

r A T “l 

g x (s) 

(4.15) 

II 

A TT 

_g (s)_ 



T 



Using (4.14), the quantities we want to estimate can be expressed in terms of 
the vector t 


(4.16) 




v?(l,s) 


L 


P (-l,s) 


P n (l,s) 


T 

L_ ~ 


II 


From formulas (4.6) and (4.16) we deduce 
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= E _1 (s) 

g z (s) 

(4.17) 




* TT 

Zn < 1>s) 


A II, \ 
& (s) 


and, since by (4.7b) E(s) has a uniformly bounded inverse, we end up with: 

< 4 * i8 > ig(-i..)i 2 + ii?<i..)i 2 < K ijg<s)i 2 , 

Re s = n > Hq > 0. 

Thus, (4.13) and (4.18) add up to (1.14), which is the definition of 
stability. 

We note that once factorized as in (4.6), the boundary matrix D(s) is 
invertible, if and only if the matrix E(s) is. Both Assumptions I and III 
were introduced in order to further guarantee that E(s) is uniformly 
bounded; in fact, one could assume, instead, that an s— uniform estimate like 
(4.7b) holds. This corresponds to the Uniform Krelss Condition which 
characterizes the stability of difference approximations to initial-boundary 
value systems. The merit of Assumption III, however, is that it deals with 
the scalar problem only, rather than with the intricate coupled boundary 
matrix E(s). 

5. THE CHEBYSHEV, LEGENDRE, AND GEGENBAUER PSEUDOSPECTRAL METHODS 

In this section we address ourselves to the question of the stability of 
the most commonly used pseudospectral methods, namely Legendre and Chebyshev 
collocation. We employ the stability criterion stated in Theorem 4.3, i.e., 
we shall verify that the scalar stability Assumptions II and III hold. 
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The scalar model is, of course: 


(5.1) 


SV N = a dF' + q(x) t(s) 


with v given at x = 1 ( — 1 ) if a > 0 (a < 0), respectively. 


Our results may be stated in a more general context, involving Gegenbauer 
polynomials, of which Chebyshev and Legendre polynomials are particular 
cases. The Gegenbauer polynomials C^ + ^ are (suitably normalized) 
orthogonal polynomials with respect to the weight 


(5.2) 


u>(x) = (1 - X 2 )“ 


in the interval -1 < x < 1. Thus, the Legendre polynomial Pjj is a multiple 

l/> 0 

of L , and the Chebyshev polynomial Tjg is a multiple of C^. 

We shall make extensive use of the Lobatto quadrature rule, which is valid 

whenever p(x) is a polynomial of degree < 2N+2: 


(5.3) 


1 N+l 

J to(x) p(x) dx = l to. p(x . ) , 

_ i a — n J J 


to . > 0. 
J 


Here, Xq = 1 and x n+i = ~1 denote the endpoints of the interval of 

integration and the interior points 


(5.4) 


1 > x^ > x 2 > 


X N— 1 > > _1 


a+ | 

are the zeros of the Gegenbauer polynomial . In particular, for the 
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Chebyshev method, these points are the zeros of , while for the Legendre 

method they are the zeros of P' . - this follows from the identity: 


(5.5) 


d r a 
dx N 


2a C 


a+1 
N-l * 


We consider collocating the equations 

x j , 1 < j < N. Then the forcing polynomial q(x) 
and may be normalized as 


(5.6) 


q(x) 






> 


at the interior points 

a+ 

is proportional to 


in order to have q(l) = 1. 


THEOREM 5.1. (Scalar stability). 

Let -1 < a < 0; then we have for all s with Re s = i) > 0: 

(5.7a) q / (1+x) io(x) |v(x,s)| 2 dx < const A(N) |v(l,s)| 2 , 

-1 

(5.7b) q / (1-x) u>(x) |v(x,s)| 2 dx < const A(N) |v(-l,s)| 2 . 


-1 

Here , A(N) is a power growth bound: 


a > 0, 


a < 0. 


/ N~ 2ct 
A(N) = j 

(n 2 


-1 < a < 0, 


(5.7c) 


a = 0 


(Legendre), 


on | cm 
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and the corresponding, weighting pair is given by 


(5.7d) 


to 1 = (1-x) (l-x 2 ) a 
= (1+x) (1-X 2 ) a . 


Proof. Let us study first the case a > 0 - without loss of generality 

we may take a = 1. Consider the Nth degree polynomial 

( 5 * 8 ) z(x,s) = v(x,s) - q(x) v(l,s). 

We have: 

A 

(5.9a) z(Xj,s) = v(Xj , s) , 1 < j < N, 

(5.9b) z(l,s) = 0, 

and, in view of (5.1): 

Jj A 

(5.9c) sz(x ,s) = — z(x ,s) + v(l ,s) q'(x. ) for 1 < j < N. 

J ax J J 

Multiply (5.9c) by (l+x^)z (x^ ,s)ojj and sum (terms with j = 0 and j = N + 1 
can be included, as their contribution vanishes) 
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(5.10) 


SUj(1+Xj) |z(Xj,s)| Z = ^ Uj(1+Xj)(^ |z| 2 )(x j>S ) 


A ^ 

+ v(l,s) l m (1+x ) z (x ,s) q'(x. ). 
j=0 J J J J 


On the left of this equality we have z(xj,s) = v(xj,s), while the sums on the 
right are in fact exact with the corresponding integrals in (5.3) and we 
obtain: 


(5.11) 

N 1 ^1 

s l (1+x. )uj. |v(x ,s| 2 = / (l+s)w — | z | 2 dx + v(l,s) • / (l+x)u)Z q' dx. 

j = l J J J _i -1 

Integrating by parts the first term on the right-hand side and using the 
Cauchy-Schwartz inequality on the second, we conclude: 

(5.12) 

N 9 

n l (1+Xj)l0j | v(Xj > s ) | = 

= (l+x)io | z j 2 - / ((l+x)m)'' | z j 2 dx + Re(v(l,s) J (l+x)mz q' dx) 

-1 -1 -1 

11 12 
< - / ((l+x)w)' I Z 1 2 + J ( ( l+x)to ) * | Z 1 2 + | v( 1 » s ) | 2 J dx 

= | v( 1 , s ) | 2 A(N), 


where 
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(5.13) 


A(N) = / 
-1 


1 ((l+x)to(x) q'(x))‘ 
[ (l+x)m(x ) )' 


dx. 


Note that by adding the boundary terms to the sum on the left, it can be 
replaced by the corresponding integral in (5.7), and that all the integrals 
involved are proper integrals, because of the limitations on a. It remains 
to show that A(N) is bounded by some power of N. A coarse bound can be 
obtained immediately, by remarking that 


(5.14) ^ max ^ |C*(x)| = C*(l), for X > 0, 

and therefore |q(x)| is bounded by 1 in -1 < x < 1. By a well known 
extremal property of Chebyshev polynomial, [6 .Theorem 2.24], the derivatives 
cannot grow faster than the corresponding derivatives of T N , and, thus, 

|q'| < N 2 . Hence, 


(5.15) 


A(N) < N F(a) , 


F(a) = / 

-1 


(l+x) a+2 (l-x) a+1 
1 + (2a-l)x 


dx* 


The more delicate estimates of (5.7c) are relegated to Appendix B. The proof 
of (5.7b) follows along similar lines. 

We now show that the collocation methods we have described satisfy 
Assumption III. 


THEOREM 5.2. Let v be a polynomial in x which satisfies 


A 1 A 

sv(x,s) = ^ v( x , s ) + q(x) 


(5.16) 
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where q(x) is defined by (5.6). If v(x=1,Sq) = 0, then Re sg < 0. _In 

A 

fact , v satisfies (4.1b), i.e.. Assumption III holds . 

Proof « In a manner similar to (5.10) - (5.12) we can deduce 

(5.17) Re s^M-x^ |v( Xj ,s 0 )| 2 = (l+x^io.. — |v| 2 (x^Sg), 0 < j < N+l , 

(5.18) 

Re s Q / (l+x)m |v| 2 dx = / (l+x)ui 4 — |v 2 | dx = - J | v | 2 ((l+x)u)' dx. 

-1 -1 -1 

From the last formula, it is obvious that Re sg < 0. Because of the symmetry 
of q(x) , 

(5.19) q(x) = q N (x) = (~1) N q N (-x), 

it is clear that v(-l,s) = (— 1)^ v(l,-s). Consider now the rational 
function of s, with real coefficients: 

A 

(5.20) f ( s ) — y. ( . 1 ».~, s > . 

v(l ,s) 

The limit of |f(s)| exists uniformly as s » and equals one, while on the 
imaginary axis 


(5.21) 


f(s) | = 1, 
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being the quotient of a number and its conjugate. Therefore, since 

a 

N v (— 1 s ) 

f(s) = (-1) ^ - — is regular for Re s > 0, it is bounded by one in 

v(l,s) 

magnitude - which is Assumption III. 

With minor modifications the same proof covers the case a < 0* 

Finally, we can summarize this section by applying Theorem 4.3 to the 
scalar results we have obtained that far. 

THEOREM 5.3. (Stability of Gegenbauer Collocation for Systems). 

Consider the hyperbolic system (1.1), satisfying Assumption I. Discretize 

it by collocating at the zeros of C N , where -1 < a < 0. The resulting 
spectral method is stable, in the sense of Lemma 1.2, with parameters defined 
by (5.7c) and (5.7d). In particular, Chebyshev and Legendre collocation 
methods are is stable. 
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APPENDIX A 

In this section we prove several a priori bounds for the solution of the 
differential system (1.1) and its scalar counterparts. 

Introducing into (1.1) the new variable e nt u(x,t) and taking the 
Fourier transform with respect to time, we obtain the equations 


(A. la) 


(A. lb) 


e AI A 1 

s u = A 


du 

A/ 

dx 


"II A II 
s u = A 


4 11 

dx 


-1 < x < 1 , 


with boundary conditions 


(A. 2) 


u I (-l,s) = L u II (-l,s) + g I (s) 
u n (l,s) = R u^ ( 1 , s ) + g n (s). 


Set a) (x) 1 + ex with 0 < e < 1 yet to be determined, and multiply 


1 -I* 


(A. la) by to u . After integration by parts, we find: 


n IIU I (x,s)ll^ I = 0) I (x) U 1 (XjS)A 1 U^XjS) 


1 f 1 AI * A 1 I J 

“ e I u A u dx, 
x=-l -1 


where p - Re s. Since eA 1 is negative, the second term on the right-hand 
side does not exceed ^ z _ | A | Hu I (x,s)l 2 I and hence 




(A.. 3) (p - p ) ||u i (x,s)ll‘ i < (l+ex)u* (x,s) A* u*(x,s) 

VJ ^ ~ r»j 

a) 


1 

x=-l 
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With *0 = 1' £ l A l* 


We can draw two conclusions from (A. 3). First, by applying it to each of the 
scalar inflow equations separately, we have for all s with Re s = n > r^, a 
positive left-hand side 


0 < a j[ (1+e) |u (j) (l,s)| 2 - (1-e) |u^(-l ,s) | 2 ] , 1 < j < l, 

that is, as aj < 0 


(A.4) 


U I (l,s) 


U^C-l.s) 


1 2 e 

1 + e 


u 1 = 


u^(x,s) 


u^ ^(x,s) 


Secondly, since (1 + ex)A 1 is negative, the inequality (A. 3) yields: 


(A. 5) (n - tIq) II u 1 Cx, s ) II 2 j < (1 - e) [A 1 1 |u^(-l,s)| 2 . 

w 


We can treat similarly the outflow part of the system. Choosing 
co 11 = 1 - ex, we find that 


(A. 6) 


U II (-1 >s) 


U n (l,s) 


1-e 

1 + e * 


,11 


*(£+!)✓ V 

U (x,s) 


U^(x,s) I 


and 
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(A.7) 


(n “ n n ) llu II (x,s)n < ( l-e) I A 11 1 |u II (l,s)| 


Finally, we rewrite the boundary conditions in the form 


i-R 


U I (-l,s) 


-L 


U I 1 (~l,s) 

U II (l,s) 



u (~l,s) 


''II/, V 

u (l,s) 


' 'T-r \ 
g (s) 


'II, ^ 


, Re s = ri > 


As in Lemma 4, the matrix on the left is invertible; in view of (a. 4) and 
(A. 6) the inverse is bounded in the norm 


-R 


U I (l , s) 
U I (-1 ,s) 


-L 


U n (-l,s) 

U n (l,s) 


-1 


1 + C|R| + |L|) f-=-| 

i - l R l N (iJf) 2 


provided e is chosen so that the denominator is positive. Then we have: 


(A. 8) |u I (-l,s)|^+ |u J " L (l,s)| z < const[ |g 1 (s) | z + Ig^Cs) | z ], 


1 2 . m 




and together with (A. 5) and (A. 7) we reach 

(A. 9) (n - rig) llu(x,s)|| 2 < const |g(s)| 2 , Re s = n > rig. 

Hence, the energy estimate referred to in (1.2b) follows from Parseval's 
relation, with any Tig satisfying 
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(A. 10) n Q - | A | > H - R -l- l . L . l) V2 ~ 1 | a | . 

2( | R| | L | ) 1/2 

We conclude by noting that in the dissipative case, where |r| |l| < 1, 
one may take n Q = e = 0 and the formulas (A. 4), (A. 6), then lead to 


u^(l,s)| < |u^\-l,s)| , 

l < j < a 

u (j) (-l,s)| < |u (j) (l,s)|, 

i < j < n, 


for Re s = n > 0. 
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APPENDIX B 

Here we evaluate (5.13), in order to obtain the asymptotic behavior of our 
estimates for large N. We shall compute 


(B.l) A(N) = / 
-1 


1 ((l+x)ujq ') 2 l [(l+x)(l-x 2 ) a (c“ + * ^ (x ))"] 2 


((l+x)io)' 


d x = J 
-1 


[a+x)u- x 2 )°]' [<c^ 1 1,2 )"(d ] 2 


dx. 


separately for a = 0, the Legendre spectral method, and for -1 < a < 0. In 
the first case, we integrate-by-parts 


(B.2) 


1 

/ (l+x) z (P^) z dx 


P N+l (x) P N+l (x)(1+X 


-1 


- I P N«( P S + 1 (I+X) )' <**• 


-1 


Once we obtain an integral involving P' + ^ > we can use Lobatto quadrature to 
exploit the fact that P' + ^ vanishes at the interior nodes. Then use the 
explicit estimates 


( 7 ^ p n) (1) ■ °<» 2 "> 

dx 

and 

u) Q = 0(N _2 ) 


to deduce that A(N) = 0(N 2 ). 


In the second case we proceed by first estimating 
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(B.3) A(N) = / 


-1 


1 (l-x) a+2 (l-x) a+1 
1 + (2a-l)x 




Kit ' 2 )'<» 


dx 


rain( |a | \\ - |«|j _{ (1+x)t<+2 (1 - x) “ 


(C ‘ /2 )"<*> 




dx. 


We then substitute from the Gegenbauer differential equation 


(B.4) (l-x 2 )(^)» ■ (2a+l) x (C^/2 )' - (K+lXMaW)^ 


to reach 


(B.5) A(N) 


(<&‘ /2 )'a> 


/ U-x 2 )“(l+x)(2a+l)x(C^ 1 1/ 2 )" dx 


-1 


- (N+l ) (N+2a+2 ) / (l-x 2 ) a (l+x) (cjj^ )» dx 


In this expression, the second integral vanishes by orthogonality, and the 
first may be evaluated by Lobatto quadrature 


(B.6) 


A(N) 


«£i 1 / 2 ra> 


2 v (C^rcDtc^ra). 


N+l 


As 


a+ 

(l+x)C N 


3 

7 


vanishes at the nodes x j , 1 < j < N+l , we can compute 
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a+ | 

As (l+x)C^ vanishes at the nodes x j , 1 < j < N+l , we can compute 


(B.7) 


“0 = 


J 

-1 


. 3 

? a+ T 

1 (1-x )° (1+x) C N (x)dx 

2c n a) 


which we estimate by the means of (5.14). Finally, we obtain 


(B.8) 


A(N) < 0(N 2ot ) 


in accordance to formula (5.7c). 
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